home *** CD-ROM | disk | FTP | other *** search
/ Magnum One / Magnum One (Mid-American Digital) (Disc Manufacturing).iso / d11 / frasrc14.arc / FPU087.ASM < prev    next >
Assembly Source File  |  1990-08-02  |  31KB  |  1,362 lines
  1. TITLE fpu087.asm (C) 1989, Mark C. Peterson, CompuServe [70441,3353]
  2. SUBTTL All rights reserved.
  3. ;
  4. ;  Code may be used in any program provided the author is credited
  5. ;    either during program execution or in the documentation.  Source
  6. ;    code may be distributed only in combination with public domain or
  7. ;    shareware source code.  Source code may be modified provided the
  8. ;    copyright notice and this message is left unchanged and all
  9. ;    modifications are clearly documented.
  10. ;
  11. ;    I would appreciate a copy of any work which incorporates this code.
  12. ;
  13. ;    Mark C. Peterson
  14. ;    253 West St., H
  15. ;    Plantsville, CT 06579
  16. ;    (203) 276-9474
  17. ;
  18. ;  References:
  19. ;     The VNR Concise Encyclopedia of Mathematics
  20. ;        by W. Gellert, H. Hustner, M. Hellwich, and H. Kastner
  21. ;        Published by Van Nostrand Reinhold Comp, 1975
  22. ;
  23. ;     80386/80286 Assembly Language Programming
  24. ;        by William H. Murray, III and Chris H. Pappas
  25. ;        Published by Osborne McGraw-Hill, 1986
  26. ;        
  27. ;
  28. ;
  29.  
  30. IFDEF ??version
  31. MASM51
  32. QUIRKS
  33. EMUL
  34. ENDIF
  35.  
  36. .model medium, c
  37.  
  38. extrn cos:far
  39. extrn _Loaded387sincos:far
  40. extrn compiled_by_turboc:word
  41.  
  42.  
  43. .data
  44.  
  45. extrn cpu:WORD
  46.  
  47. PUBLIC TrigLimit, TrigOverflow
  48.  
  49. PiFg13         dw       6487h
  50. InvPiFg33      dd       0a2f9836eh
  51. InvPiFg16      dw       517ch
  52. Ln2Fg16        dw       0b172h
  53. TrigOverflow   dw       0
  54. TrigLimit      dd       0
  55. one            dw       ?
  56. expSign        dw       ?
  57. a              dw       ?
  58. SinNeg         dw       ?
  59. CosNeg           dw    ?
  60. Ans           dq    ?
  61.  
  62.  
  63. TaylorTerm  MACRO
  64. LOCAL Ratio
  65.    add   Factorial, one
  66.    jnc   SHORT Ratio
  67.  
  68.    rcr   Factorial, 1
  69.    shr   Num, 1
  70.    shr   one, 1
  71.  
  72. Ratio:
  73.    mul   Num
  74.    div   Factorial
  75. ENDM
  76.  
  77.  
  78.  
  79. _4_         dq    4.0
  80. _2_         dq    2.0
  81. _1_         dq    1.0
  82. PointFive   dq    0.5
  83. temp        dq     ?
  84. Sign        dw     ?
  85.  
  86. extrn fpu:word
  87.  
  88. .code
  89.  
  90.  
  91. FPUcplxmul     PROC     x:word, y:word, z:word
  92.    mov   bx, x
  93.    fld   QWORD PTR [bx]       ; x.x
  94.    fld   QWORD PTR [bx+8]     ; x.y, x.x
  95.    mov   bx, y
  96.    fld   QWORD PTR [bx]       ; y.x, x.y, x.x
  97.    fld   QWORD PTR [bx+8]     ; y.y, y.x, x.y, x.x
  98.    mov   bx, z
  99.    fld   st                   ; y.y, y.y, y.x, x.y, x.x
  100.    fmul  st, st(3)            ; y.y*x.y, y.y. y.x, x.y, x.x
  101.    fld   st(2)                ; y.x, y.y*x.y, y.y, y.x, x.y, x.x
  102.    fmul  st, st(5)            ; y.x*x.x, y.y*x.y, y.y, y.x, x.y, x.x
  103.    fsubr                      ; y.x*x.x - y.y*x.y, y.y, y.x, x.y, x.x
  104.    fstp  QWORD PTR [bx]       ; y.y, y.x, x.y, x.x
  105.    fmulp st(3), st            ; y.x, x.y, x.x*y.y
  106.    fmul                       ; y.x*x.y, x.x*y.y
  107.    fadd                       ; y.x*x.y + x.x*y.y
  108.    fstp  QWORD PTR [bx+8]
  109.    ret
  110. FPUcplxmul     ENDP
  111.  
  112.  
  113.  
  114.  
  115.  
  116. FPUcplxdiv     PROC     x:word, y:word, z:word
  117.    mov   bx, x
  118.    fld   QWORD PTR [bx]       ; x.x
  119.    fld   QWORD PTR [bx+8]     ; x.y, x.x
  120.    mov   bx, y
  121.    fld   QWORD PTR [bx]       ; y.x, x.y, x.x
  122.    fld   QWORD PTR [bx+8]     ; y.y, y.x, x.y, x.x
  123.    fld   st                   ; y.y, y.y, y.x, x.y, x.x
  124.    fmul  st, st               ; y.y*y.y, y.y, y.x, x.y, x.x
  125.    fld   st(2)                ; y.x, y.y*y.y, y.y, y.x, x.y, x.x
  126.    fmul  st, st               ; y.x*y.x, y.y*y.y, y.y, y.x, x.y, x.x
  127.    fadd                       ; mod, y.y, y.x, x.y, x.x
  128.    fdiv  st(1), st            ; mod, y.y=y.y/mod, y.x, x.y, x.x
  129.    fdivp st(2), st            ; y.y, y.x=y.x/mod, x.y, x.x
  130.    mov   bx, z
  131.    fld   st                   ; y.y, y.y, y.x, x.y, x.x
  132.    fmul  st, st(3)            ; y.y*x.y, y.y. y.x, x.y, x.x
  133.    fld   st(2)                ; y.x, y.y*x.y, y.y, y.x, x.y, x.x
  134.    fmul  st, st(5)            ; y.x*x.x, y.y*x.y, y.y, y.x, x.y, x.x
  135.    fadd                       ; y.x*x.x - y.y*x.y, y.y, y.x, x.y, x.x
  136.    fstp  QWORD PTR [bx]       ; y.y, y.x, x.y, x.x
  137.    fmulp st(3), st            ; y.x, x.y, x.x*y.y
  138.    fmul                       ; y.x*x.y, x.x*y.y
  139.    fsubr                      ; y.x*x.y + x.x*y.y
  140.    fstp  QWORD PTR [bx+8]
  141.    ret
  142. FPUcplxdiv     ENDP
  143.  
  144.  
  145.  
  146. FPUcplxlog     PROC     x:word, z:word
  147. LOCAL Status:word
  148.    mov   bx, x
  149.    fld   QWORD PTR [bx+8]        ; x.y
  150.    fld   QWORD PTR [bx]          ; x.x, x.y
  151.    mov   bx, z
  152.    fldln2                        ; ln2, x.x, x.y
  153.    fdiv  _2_                     ; ln2/2, x.x, x.y
  154.    fld   st(2)                   ; x.y, ln2/2, x.x, x.y
  155.    fmul  st, st                  ; sqr(x.y), ln2/2, x.x, x.y
  156.    fld   st(2)                   ; x.x, sqr(x.y), ln2/2, x.x, x.y
  157.    fmul  st, st                  ; sqr(x.x), sqr(x.y), ln2/2, x.x, x.y
  158.    fadd                          ; mod, ln2/2, x.x, x.y
  159.    fyl2x                         ; z.x, x.x, x.y
  160.    fstp  QWORD PTR [bx]          ; x.x, x.y
  161.    cmp   fpu, 387
  162.    jne   Restricted
  163.  
  164.    fpatan
  165.    jmp   StoreZX
  166.  
  167. Restricted:
  168.    mov   bx, x
  169.    mov   dh, BYTE PTR [bx+7]
  170.    or    dh, dh
  171.    jns   ChkYSign
  172.  
  173.    fchs                          ; |x.x|, x.y
  174.  
  175. ChkYSign:
  176.    mov   dl, BYTE PTR [bx+8+7]
  177.    or    dl, dl
  178.    jns   ChkMagnitudes
  179.  
  180.    fxch                          ; x.y, |x.x|
  181.    fchs                          ; |x.y|, |x.x|
  182.    fxch                          ; |x.x|, |x.y|
  183.  
  184. ChkMagnitudes:
  185.    fcom  st(1)                   ; x.x, x.y
  186.    fstsw Status                  ; x.x, x.y
  187.    test  Status, 4500h
  188.    jz    XisGTY
  189.  
  190.    test  Status, 4000h
  191.    jz    XneY
  192.  
  193.    fstp  st                      ; x.y
  194.    fstp  st                      ; <empty>
  195.    fldpi                         ; Pi
  196.    fdiv  _4_                     ; Pi/4
  197.    jmp   ChkSignZ
  198.  
  199. XneY:
  200.    fxch                          ; x.y, x.x
  201.    fpatan                        ; Pi/2 - Angle
  202.    fldpi                         ; Pi, Pi/2 - Angle
  203.    fdiv  _2_                     ; Pi/2, Pi/2 - Angle
  204.    fsubr                         ; Angle
  205.    jmp   ChkSignZ
  206.  
  207. XisGTY:
  208.    fpatan
  209.  
  210. ChkSignZ:
  211.    or    dh, dh
  212.    js    NegX
  213.  
  214.    or    dl, dl
  215.    jns   StoreZX
  216.  
  217.    fchs
  218.    jmp   StoreZX
  219.  
  220. NegX:
  221.    or    dl, dl
  222.    js    QuadIII
  223.  
  224.    fldpi
  225.    fsubr
  226.    jmp   StoreZX
  227.  
  228. QuadIII:
  229.    fldpi
  230.    fsubr
  231.    fchs
  232.  
  233. StoreZX:
  234.    mov   bx, z
  235.    fstp  QWORD PTR [bx+8]        ; <empty>
  236.    ret
  237. FPUcplxlog     ENDP
  238.  
  239.  
  240.  
  241.  
  242. FPUsinhcosh    PROC     x:word, sinh:word, cosh:word
  243. LOCAL Control:word
  244.    fstcw Control
  245.    push  Control                       ; Save control word on the stack
  246.    or    Control, 0000110000000000b 
  247.    fldcw Control                       ; Set control to round towards zero
  248.  
  249.    mov   Sign, 0              ; Assume the sign is positive
  250.    mov   bx, x
  251.  
  252.    fldln2                     ; ln(2)
  253.    fdivr QWORD PTR [bx]       ; x/ln(2)
  254.  
  255.    cmp   BYTE PTR [bx+7], 0
  256.    jns   DuplicateX
  257.  
  258.    fchs                       ; x = |x|
  259.  
  260. DuplicateX:   
  261.    fld   st                   ; x/ln(2), x/ln(2)
  262.    frndint                    ; int = integer(|x|/ln(2)), x/ln(2)
  263.    fxch                       ; x/ln(2), int
  264.    fsub  st, st(1)            ; rem < 1.0, int
  265.    fdiv  _2_                  ; rem/2 < 0.5, int
  266.    f2xm1                      ; (2**rem/2)-1, int
  267.    fadd  _1_                  ; 2**rem/2, int
  268.    fmul  st, st               ; 2**rem, int
  269.    fscale                     ; e**|x|, int
  270.    fstp  st(1)                ; e**|x|
  271.  
  272.    cmp   BYTE PTR [bx+7], 0
  273.    jns   ExitFexp
  274.  
  275.    fdivr _1_                  ; e**x      
  276.  
  277. ExitFexp:
  278.    fld   st                   ; e**x, e**x
  279.    fdivr PointFive            ; e**-x/2, e**x
  280.    fld   st                   ; e**-x/2, e**-x/2, e**x
  281.    fxch  st(2)                ; e**x, e**-x/2, e**-x/2
  282.    fdiv  _2_                  ; e**x/2,  e**-x/2, e**-x/2
  283.    fadd  st(2), st            ; e**x/2,  e**-x/2, cosh(x)
  284.    fsubr                      ; sinh(x), cosh(x)
  285.  
  286.    mov   bx, sinh             ; sinh, cosh
  287.    fstp  QWORD PTR [bx]       ; cosh
  288.    mov   bx, cosh
  289.    fstp  QWORD PTR [bx]       ; <empty>
  290.  
  291.    pop   Control
  292.    fldcw Control              ; Restore control word
  293.    ret
  294. FPUsinhcosh    ENDP
  295.  
  296.  
  297. FPUsincos  PROC  x:word, sinx:word, cosx:word
  298. LOCAL Status:word
  299.    mov   bx, x
  300.    fld   QWORD PTR [bx]       ; x
  301.  
  302.    cmp   fpu, 387
  303.    jne   Use387FPUsincos
  304.  
  305.    call  _Loaded387sincos     ; cos(x), sin(x)
  306.    mov   bx, cosx
  307.    fstp  QWORD PTR [bx]       ; sin(x)
  308.    mov   bx, sinx
  309.    fstp  QWORD PTR [bx]       ; <empty>
  310.    ret
  311.  
  312. Use387FPUsincos:
  313.  
  314.    sub   sp, 8                ; save 'x' on the CPU stack
  315.    mov   bx, sp
  316.    fstp  QWORD PTR [bx]       ; FPU stack:  <empty>
  317.  
  318.    call  cos
  319.  
  320.    add   sp, 8                ; take 'cos(x)' off the CPU stack
  321.    mov   bx, ax
  322.    cmp   compiled_by_turboc,0
  323.    jne   turbo_c1
  324.  
  325.    fld   QWORD PTR [bx]       ; FPU stack:  cos(x)
  326.  
  327. turbo_c1:
  328.    fld   st                   ; FPU stack:  cos(x), cos(x)
  329.    fmul  st, st               ; cos(x)**2, cos(x)
  330.    fsubr _1_                  ; sin(x)**2, cos(x)
  331.    fsqrt                      ; +/-sin(x), cos(x)
  332.  
  333.    mov   bx, x
  334.    fld   QWORD PTR [bx]       ; x, +/-sin(x), cos(x)
  335.    fldpi                      ; Pi, x, +/-sin(x), cos(x)
  336.    fadd  st, st               ; 2Pi, x, +/-sin(x), cos(x)
  337.    fxch                       ; |x|, 2Pi, +/-sin(x), cos(x)
  338.    fprem                      ; Angle, 2Pi, +/-sin(x), cos(x)
  339.    fstp  st(1)                ; Angle, +/-sin(x), cos(x)
  340.    fldpi                      ; Pi, Angle, +/-sin(x), cos(x)
  341.  
  342.    cmp   BYTE PTR [bx+7], 0
  343.    jns   SignAlignedPi
  344.  
  345.    fchs                       ; -Pi, Angle, +/-sin(x), cos(x)
  346.  
  347. SignAlignedPi:
  348.    fcompp                     ; +/-sin(x), cos(x)
  349.    fstsw Status               ; +/-sin(x), cos(x)
  350.  
  351.    mov   ax, Status
  352.    and   ah, 1
  353.    jz    StoreSinCos          ; Angle <= Pi
  354.  
  355.    fchs                       ; sin(x), cos(x)
  356.  
  357. StoreSinCos:
  358.    mov   bx, sinx
  359.    fstp  QWORD PTR [bx]       ; cos(x)
  360.    mov   bx, cosx
  361.    fstp  QWORD PTR [bx]       ; <empty>
  362.    ret
  363. FPUsincos   ENDP
  364.  
  365.  
  366. PUBLIC r16Mul
  367. r16Mul     PROC    uses si di, x1:word, x2:word, y1:word, y2:word
  368.       mov   si, x1
  369.       mov   bx, x2
  370.       mov   di, y1
  371.       mov   cx, y2
  372.  
  373.       xor   ax, ax
  374.       shl   bx, 1
  375.       jz    Exitr16Mult          ; Destination is zero
  376.  
  377.       rcl   ah, 1
  378.       shl   cx, 1
  379.       jnz   Chkr16Exp
  380.       xor   bx, bx               ; Source is zero
  381.       xor   si, si
  382.       jmp   Exitr16Mult
  383.  
  384.    Chkr16Exp:
  385.       rcl   al, 1
  386.       xor   ah, al               ; Resulting sign in ah
  387.       stc                        ; Put 'one' bit back into number
  388.       rcr   bl, 1
  389.       stc
  390.       rcr   cl, 1
  391.  
  392.       sub   ch, 127              ; Determine resulting exponent
  393.       add   bh, ch
  394.       mov   al, bh
  395.       mov   es, ax               ; es has the resulting exponent and sign
  396.  
  397.       mov   ax, di
  398.       mov   al, ah
  399.       mov   ah, cl
  400.  
  401.       mov   dx, si
  402.       mov   dl, dh
  403.       mov   dh, bl
  404.  
  405.       mul   dx
  406.       mov   cx, es
  407.  
  408.       shl   ax, 1
  409.       rcl   dx, 1
  410.       jnc   Remr16MulOneBit      ; 'One' bit is the next bit over
  411.  
  412.       inc   cl                   ; 'One' bit removed with previous shift
  413.       jmp   Afterr16MulNorm
  414.  
  415.    Remr16MulOneBit:
  416.       shl   ax, 1
  417.       rcl   dx, 1
  418.  
  419.    Afterr16MulNorm:
  420.       mov   bl, dh               ; Perform remaining 8 bit shift
  421.       mov   dh, dl
  422.       mov   dl, ah
  423.       mov   si, dx
  424.       mov   bh, cl               ; Put in the exponent
  425.       rcr   ch, 1                ; Get the sign
  426.       rcr   bx, 1                ; Normalize the result
  427.       rcr   si, 1
  428.    Exitr16Mult:
  429.       mov   ax, si
  430.       mov   dx, bx
  431.       ret
  432. r16Mul      ENDP
  433.  
  434.  
  435. PUBLIC RegFloat2Fg
  436. RegFloat2Fg     PROC    x1:word, x2:word, Fudge:word
  437.       mov   ax, WORD PTR x1
  438.       mov   dx, WORD PTR x2
  439.       mov   bx, ax
  440.       or    bx, dx
  441.       jz    ExitRegFloat2Fg
  442.  
  443.       xor   bx, bx
  444.       mov   cx, bx
  445.  
  446.       shl   ax, 1
  447.       rcl   dx, 1
  448.       rcl   bx, 1                   ; bx contains the sign
  449.  
  450.       xchg  cl, dh                  ; cx contains the exponent
  451.  
  452.       stc                           ; Put in the One bit
  453.       rcr   dl, 1
  454.       rcr   ax, 1
  455.  
  456.       sub   cx, 127 + 23
  457.       add   cx, Fudge
  458.       jz    ChkFgSign
  459.       jns   ShiftFgLeft
  460.  
  461.       neg   cx
  462.    ShiftFgRight:
  463.       shr   dx, 1
  464.       rcr   ax, 1
  465.       loop  ShiftFgRight
  466.       jmp   ChkFgSign
  467.  
  468.    ShiftFgLeft:
  469.       shl   ax, 1
  470.       rcl   dx, 1
  471.       loop  ShiftFgLeft
  472.  
  473.    ChkFgSign:
  474.       or    bx, bx
  475.       jz    ExitRegFloat2Fg
  476.  
  477.       not   ax
  478.       not   dx
  479.       add   ax, 1
  480.       adc   dx, 0
  481.  
  482.    ExitRegFloat2Fg:
  483.       ret
  484. RegFloat2Fg    ENDP
  485.  
  486.  
  487.  
  488. PUBLIC RegFg2Float
  489. RegFg2Float     PROC   x1:word, x2:word, FudgeFact:byte
  490.       mov   ax, x1
  491.       mov   dx, x2
  492.  
  493.       mov   cx, ax
  494.       or    cx, dx
  495.       jz    ExitFudgedToRegFloat
  496.  
  497.       mov   ch, 127 + 32
  498.       sub   ch, FudgeFact
  499.       xor   cl, cl
  500.       shl   ax, 1       ; Get the sign bit
  501.       rcl   dx, 1
  502.       jnc   FindOneBit
  503.  
  504.       inc   cl          ; Fudged < 0, convert to postive
  505.       not   ax
  506.       not   dx
  507.       add   ax, 1
  508.       adc   dx, 0
  509.  
  510.    FindOneBit:
  511.       shl   ax, 1
  512.       rcl   dx, 1
  513.       dec   ch
  514.       jnc   FindOneBit
  515.       dec   ch
  516.  
  517.       mov   al, ah
  518.       mov   ah, dl
  519.       mov   dl, dh
  520.       mov   dh, ch
  521.  
  522.       shr   cl, 1       ; Put sign bit in
  523.       rcr   dx, 1
  524.       rcr   ax, 1
  525.  
  526.    ExitFudgedToRegFloat:
  527.       ret
  528. RegFg2Float      ENDP
  529.  
  530.  
  531. PUBLIC ExpFudged
  532. ExpFudged      PROC    uses si, x_low:word, x_high:word, Fudge:word
  533. LOCAL exp:WORD
  534.       xor   ax, ax
  535.       mov   WORD PTR Ans, ax
  536.       mov   WORD PTR Ans + 2, ax
  537.       mov   ax, WORD PTR x_low
  538.       mov   dx, WORD PTR x_high
  539.       or    dx, dx
  540.       js    NegativeExp
  541.  
  542.       div   Ln2Fg16
  543.       mov   exp, ax
  544.       or    dx, dx
  545.       jz    Raiseexp
  546.  
  547.       mov   ax, dx
  548.       mov   si, dx
  549.       mov   bx, 1
  550.  
  551.    PosExpLoop:
  552.       add   WORD PTR Ans, ax
  553.       adc   WORD PTR Ans + 2, 0
  554.       inc   bx
  555.       mul   si
  556.       mov   ax, dx
  557.       xor   dx, dx
  558.       div   bx
  559.       or    ax, ax
  560.       jnz   PosExpLoop
  561.  
  562.    Raiseexp:
  563.       inc   WORD PTR Ans + 2
  564.       mov   ax, WORD PTR Ans
  565.       mov   dx, WORD PTR Ans + 2
  566.       mov   cx, -16
  567.       add   cx, Fudge
  568.       add   cx, exp
  569.       or    cx, cx
  570.       jz    ExitExpFudged
  571.       jns   LeftShift
  572.       neg   cx
  573.  
  574.    RightShift:
  575.       shr   dx, 1
  576.       rcr   ax, 1
  577.       loop  RightShift
  578.       jmp   ExitExpFudged
  579.  
  580.    NegativeExp:
  581.       not   ax
  582.       not   dx
  583.       add   ax, 1
  584.       adc   dx, 0
  585.       div   Ln2Fg16
  586.       neg   ax
  587.       mov   exp, ax
  588.  
  589.       or    dx, dx
  590.       jz    Raiseexp
  591.  
  592.       mov   ax, dx
  593.       mov   si, dx
  594.       mov   bx, 1
  595.  
  596.    NegExpLoop:
  597.       sub   WORD PTR Ans, ax
  598.       sbb   WORD PTR Ans + 2, 0
  599.       inc   bx
  600.       mul   si
  601.       mov   ax, dx
  602.       xor   dx, dx
  603.       div   bx
  604.       or    ax, ax
  605.       jz    Raiseexp
  606.  
  607.       add   WORD PTR Ans, ax
  608.       adc   WORD PTR Ans + 2, 0
  609.       inc   bx
  610.       mul   si
  611.       mov   ax, dx
  612.       xor   dx, dx
  613.       div   bx
  614.       or    ax, ax
  615.       jnz   NegExpLoop
  616.       jmp   Raiseexp
  617.  
  618.    LeftShift:
  619.       shl   ax, 1
  620.       rcl   dx, 1
  621.       loop  LeftShift
  622.  
  623.    ExitExpFudged:
  624.       ret
  625. ExpFudged      ENDP
  626.  
  627.  
  628.  
  629. PUBLIC   LogFudged
  630. LogFudged      PROC    uses si di es, x_low:word, x_high:word, Fudge:word
  631. LOCAL exp:WORD
  632.       xor   bx, bx
  633.       mov   cx, 16
  634.       sub   cx, Fudge
  635.       mov   ax, x_low
  636.       mov   dx, x_high
  637.  
  638.       or    dx, dx
  639.       jz    ChkLowWord
  640.  
  641.    Incexp:
  642.       shr   dx, 1
  643.       jz    DetermineOper
  644.       rcr   ax, 1
  645.       inc   cx
  646.       jmp   Incexp
  647.  
  648.    ChkLowWord:
  649.       or    ax, ax
  650.       jnz   Decexp
  651.       jmp   ExitLogFudged
  652.  
  653.    Decexp:
  654.       dec   cx                      ; Determine power of two
  655.       shl   ax, 1
  656.       jnc   Decexp
  657.  
  658.    DetermineOper:
  659.       mov   exp, cx
  660.       mov   si, ax                  ; si =: x + 1
  661.       shr   si, 1
  662.       stc
  663.       rcr   si, 1
  664.       mov   dx, ax
  665.       xor   ax, ax
  666.       shr   dx, 1
  667.       rcr   ax, 1
  668.       shr   dx, 1
  669.       rcr   ax, 1                   ; dx:ax = x - 1
  670.       div   si
  671.  
  672.       mov   bx, ax                  ; ax, Fudged 16, max of 0.3333333
  673.       shl   ax, 1                   ; x = (x - 1) / (x + 1), Fudged 16
  674.       mul   ax
  675.       shl   ax, 1
  676.       rcl   dx, 1
  677.       mov   ax, dx                  ; dx:ax, Fudged 35, max = 0.1111111
  678.       mov   si, ax                  ; si = (ax * ax), Fudged 19
  679.  
  680.       mov   ax, bx
  681.    ; bx is the accumulator, First term is x
  682.       mul   si                      ; dx:ax, Fudged 35, max of 0.037037
  683.       mov   es, dx                  ; Save high word, Fudged (35 - 16) = 19
  684.       mov   di, 0c000h              ; di, 3 Fudged 14
  685.       div   di                      ; ax, Fudged (36 - 14) = 21
  686.       or    ax, ax
  687.       jz    Addexp
  688.  
  689.       mov   cl, 5
  690.       shr   ax, cl
  691.       add   bx, ax                  ; bx, max of 0.345679
  692.    ; x = x + x**3/3
  693.  
  694.       mov   ax, es                  ; ax, Fudged 19
  695.       mul   si                      ; dx:ax, Fudged 38, max of 0.004115
  696.       mov   es, dx                  ; Save high word, Fudged (38 - 16) = 22
  697.       mov   di, 0a000h              ; di, 5 Fudged 13
  698.       div   di                      ; ax, Fudged (38 - 13) = 25
  699.       or    ax, ax
  700.       jz    Addexp
  701.  
  702.       mov   cl, 9
  703.       shr   ax, cl
  704.       add   bx, ax
  705.    ; x = x + x**3/3 + x**5/5
  706.  
  707.       mov   ax, es                  ; ax, Fudged 22
  708.       mul   si                      ; dx:ax, Fudged 41, max of 0.0004572
  709.       mov   di, 0e000h              ; di, 7 Fudged 13
  710.       div   di                      ; ax, Fudged (41 - 13) = 28
  711.       mov   cl, 12
  712.       shr   ax, cl
  713.       add   bx, ax
  714.  
  715.    Addexp:
  716.       shl   bx, 1                   ; bx *= 2, Fudged 16, max of 0.693147
  717.    ; x = 2 * (x + x**3/3 + x**5/5 + x**7/7)
  718.       mov   cx, exp
  719.       mov   ax, Ln2Fg16            ; Answer += exp * Ln2Fg16
  720.       or    cx, cx
  721.       js    SubFromAns
  722.  
  723.       mul   cx
  724.       add   ax, bx
  725.       adc   dx, 0
  726.       jmp   ExitLogFudged
  727.  
  728.    SubFromAns:
  729.       neg   cx
  730.       mul   cx
  731.       xor   cx, cx
  732.       xchg  cx, dx
  733.       xchg  bx, ax
  734.       sub   ax, bx
  735.       sbb   dx, cx
  736.  
  737.    ExitLogFudged:
  738.    ; x = 2 * (x + x**3/3 + x**5/5 + x**7/7) + (exp * Ln2Fg16)
  739.       ret
  740. LogFudged      ENDP
  741.  
  742.  
  743.  
  744.  
  745. PUBLIC LogFloat14
  746. LogFloat14     PROC     x1:word, x2:word
  747.       mov   ax, WORD PTR x1
  748.       mov   dx, WORD PTR x2
  749.       shl   ax, 1
  750.       rcl   dx, 1
  751.       xor   cx, cx
  752.       xchg  cl, dh
  753.  
  754.       stc
  755.         rcr   dl, 1
  756.         rcr   ax, 1
  757.  
  758.       sub   cx, 127 + 23
  759.       neg   cx
  760.       push  cx
  761.       push  dx
  762.       push  ax
  763.       call  LogFudged
  764.       add   sp, 6
  765.       ret
  766. LogFloat14     ENDP
  767.  
  768.  
  769. PUBLIC RegSftFloat
  770. RegSftFloat     PROC   x1:word, x2:word, Shift:byte
  771.       mov   ax, x1
  772.       mov   dx, x2
  773.  
  774.       shl   dx, 1
  775.       rcl   cl, 1
  776.  
  777.       add   dh, Shift
  778.  
  779.       shr   cl, 1
  780.       rcr   dx, 1
  781.  
  782.       ret
  783. RegSftFloat      ENDP
  784.  
  785.  
  786.  
  787.  
  788. PUBLIC RegDivFloat
  789. RegDivFloat     PROC  uses si di, x1:word, x2:word, y1:word, y2:word
  790.       mov   si, x1
  791.       mov   bx, x2
  792.       mov   di, y1
  793.       mov   cx, y2
  794.  
  795.       xor   ax, ax
  796.       shl   bx, 1
  797.       jnz   ChkOtherOp
  798.       jmp   ExitRegDiv           ; Destination is zero
  799.  
  800.    ChkOtherOp:
  801.       rcl   ah, 1
  802.       shl   cx, 1
  803.       jnz   ChkDivExp
  804.       xor   bx, bx               ; Source is zero
  805.       xor   si, si
  806.       jmp   ExitRegDiv
  807.  
  808.    ChkDivExp:
  809.       rcl   al, 1
  810.       xor   ah, al               ; Resulting sign in ah
  811.       stc                        ; Put 'one' bit back into number
  812.       rcr   bl, 1
  813.       stc
  814.       rcr   cl, 1
  815.  
  816.       sub   ch, 127              ; Determine resulting exponent
  817.       sub   bh, ch
  818.       mov   al, bh
  819.       mov   es, ax               ; es has the resulting exponent and sign
  820.  
  821.       mov   ax, si               ; 8 bit shift, bx:si moved to dx:ax
  822.       mov   dh, bl
  823.       mov   dl, ah
  824.       mov   ah, al
  825.       xor   al, al
  826.  
  827.       mov   bh, cl               ; 8 bit shift, cx:di moved to bx:cx
  828.       mov   cx, di
  829.       mov   bl, ch
  830.       mov   ch, cl
  831.       xor   cl, cl
  832.  
  833.       shr   dx, 1
  834.       rcr   ax, 1
  835.  
  836.       div   bx
  837.       mov   si, dx               ; Save (and shift) remainder
  838.       mov   dx, cx               ; Save the quess
  839.       mov   cx, ax
  840.       mul   dx                   ; Mult quess times low word
  841.       xor   di, di
  842.       sub   di, ax               ; Determine remainder
  843.       sbb   si, dx
  844.       mov   ax, di
  845.       mov   dx, si
  846.       jc    RemainderNeg
  847.  
  848.       xor   di, di
  849.       jmp   GetNextDigit
  850.  
  851.    RemainderNeg:
  852.       mov   di, 1                ; Flag digit as negative
  853.       not   ax                   ; Convert remainder to positive
  854.       not   dx
  855.       add   ax, 1
  856.       adc   dx, 0
  857.  
  858.    GetNextDigit:
  859.       shr   dx, 1
  860.       rcr   ax, 1
  861.       div   bx
  862.       xor   bx, bx
  863.       shl   dx, 1
  864.       rcl   ax, 1
  865.       rcl   bl, 1                ; Save high bit
  866.  
  867.       mov   dx, cx               ; Retrieve first digit
  868.       or    di, di
  869.       jz    RemoveDivOneBit
  870.  
  871.       neg   ax                   ; Digit was negative
  872.       neg   bx
  873.       dec   dx
  874.  
  875.    RemoveDivOneBit:
  876.       add   dx, bx
  877.       mov   cx, es
  878.       shl   ax, 1
  879.       rcl   dx, 1
  880.       jc    AfterDivNorm
  881.  
  882.       dec   cl
  883.       shl   ax, 1
  884.       rcl   dx, 1
  885.  
  886.    AfterDivNorm:
  887.       mov   bl, dh               ; Perform remaining 8 bit shift
  888.       mov   dh, dl
  889.       mov   dl, ah
  890.       mov   si, dx
  891.       mov   bh, cl               ; Put in the exponent
  892.       shr   ch, 1                ; Get the sign
  893.       rcr   bx, 1                ; Normalize the result
  894.       rcr   si, 1
  895.  
  896.    ExitRegDiv:
  897.       mov   ax, si
  898.       mov   dx, bx
  899.       ret
  900. RegDivFloat      ENDP
  901.  
  902.  
  903.  
  904. Term        equ      <ax>
  905. Num         equ      <bx>
  906. Factorial   equ      <cx>
  907. Sin         equ      <si>
  908. Cos         equ      <di>
  909. e           equ      <si>
  910. Inve        equ      <di>
  911.          
  912. SinCos086   PROC     uses si di, LoNum:WORD, HiNum:WORD, SinAddr:WORD, \
  913.                                 CosAddr:WORD
  914.    mov   ax, LoNum 
  915.    mov   dx, HiNum
  916.    
  917.    xor   cx, cx
  918.    mov   SinNeg, cx
  919.    mov   CosNeg, cx
  920.    mov   a, cx
  921.    or    dx, dx
  922.    jns   AnglePositive
  923.    
  924.    not   ax
  925.    not   dx
  926.    add   ax, 1
  927.    adc   dx, cx
  928.    mov   SinNeg, 1
  929.       
  930. AnglePositive:
  931.    mov   si, ax
  932.    mov   di, dx
  933.    mul   WORD PTR InvPiFg33
  934.    mov   bx, dx
  935.    mov   ax, di
  936.    mul   WORD PTR InvPiFg33
  937.    add   bx, ax
  938.    adc   cx, dx
  939.    mov   ax, si
  940.    mul   WORD PTR InvPiFg33 + 2
  941.    add   bx, ax
  942.    adc   cx, dx
  943.    mov   ax, di
  944.    mul   WORD PTR InvPiFg33 + 2
  945.    add   ax, cx
  946.    adc   dx, 0
  947.  
  948.    and   dx, 3
  949.    mov   a, dx
  950.  
  951.    mov   Num, ax
  952.    mov   Factorial, WORD PTR InvPiFg33 + 2
  953.    mov   one, Factorial
  954.    mov   Cos, Factorial          ; Cos = 1
  955.    mov   Sin, Num                ; Sin = Num
  956.       
  957. LoopIntSinCos:
  958.    TaylorTerm                    ; Term = Num * (x/2) * (x/3) * (x/4) * . . .
  959.    sub   Cos, Term               ; Cos = 1 - Num*(x/2) + (x**4)/4! - . . .
  960.    cmp   Term, WORD PTR TrigLimit
  961.    jbe   SHORT ExitIntSinCos
  962.  
  963.    TaylorTerm
  964.    sub   Sin, Term               ; Sin = Num - Num*(x/2)*(x/3) + (x**5)/5! - . . .
  965.    cmp   Term, WORD PTR TrigLimit
  966.    jbe   SHORT ExitIntSinCos
  967.       
  968.    TaylorTerm
  969.    add   Cos, Term
  970.    cmp   Term, WORD PTR TrigLimit
  971.    jbe   SHORT ExitIntSinCos
  972.       
  973.    TaylorTerm                    ; Term = Num * (x/2) * (x/3) * . . .
  974.    add   Sin, Term
  975.    cmp   Term, WORD PTR TrigLimit
  976.    jnbe  LoopIntSinCos
  977.       
  978. ExitIntSinCos:
  979.    xor   ax, ax
  980.    mov   cx, ax
  981.    cmp   Cos, WORD PTR InvPiFg33 + 2
  982.    jb    CosDivide               ; Cos < 1.0
  983.       
  984.    inc   cx                      ; Cos == 1.0
  985.    jmp   StoreCos
  986.       
  987. CosDivide:
  988.    mov   dx, Cos
  989.    div   WORD PTR InvPiFg33 + 2
  990.       
  991. StoreCos:
  992.    mov   Cos, ax                 ; cx:Cos
  993.  
  994.    xor   ax, ax
  995.    mov   bx, ax
  996.    cmp   Sin, WORD PTR InvPiFg33 + 2
  997.    jb    SinDivide               ; Sin < 1.0
  998.    
  999.    inc   bx                      ; Sin == 1.0
  1000.    jmp   StoreSin
  1001.       
  1002. SinDivide:
  1003.    mov   dx, Sin
  1004.    div   WORD PTR InvPiFg33 + 2
  1005.       
  1006. StoreSin:
  1007.    mov   Sin, ax                 ; bx:Sin
  1008.  
  1009.    test  a, 1
  1010.    jz    ChkNegCos
  1011.  
  1012.    xchg  cx, bx
  1013.    xchg  Sin, Cos
  1014.    mov   ax, SinNeg
  1015.    xchg  ax, CosNeg
  1016.    mov   CosNeg, ax
  1017.  
  1018. ChkNegCos:
  1019.    mov   ax, a
  1020.    shr   al, 1
  1021.    rcl   ah, 1
  1022.    xor   ah, al
  1023.    jz    ChkNegSin
  1024.  
  1025.    xor   CosNeg, 1
  1026.  
  1027. ChkNegSin:
  1028.    test  a, 2
  1029.    jz    CorrectQuad
  1030.  
  1031.    xor   SinNeg, 1
  1032.  
  1033. CorrectQuad:
  1034.  
  1035.    cmp   CosNeg, 1
  1036.    jne   CosPolarized
  1037.  
  1038.    not   Cos
  1039.    not   cx 
  1040.    add   Cos, 1
  1041.    adc   cx, 0
  1042.  
  1043. CosPolarized:     
  1044.    mov   dx, bx
  1045.    mov   bx, CosAddr
  1046.    mov   WORD PTR [bx], Cos
  1047.    mov   WORD PTR [bx+2], cx
  1048.  
  1049.    cmp   SinNeg, 1
  1050.    jne   SinPolarized
  1051.  
  1052.    not   Sin
  1053.    not   dx
  1054.    add   Sin, 1
  1055.    adc   dx, 0
  1056.  
  1057. SinPolarized:
  1058.    mov   bx, SinAddr
  1059.    mov   WORD PTR [bx], Sin
  1060.    mov   WORD PTR [bx+2], dx
  1061.    ret
  1062. SinCos086      ENDP
  1063.       
  1064.       
  1065.       
  1066. _e2x   PROC
  1067.    mov   expSign, 0
  1068.    or    dx, dx
  1069.    jns   CalcExp
  1070.       
  1071.    mov   expSign, 1
  1072.    not   ax
  1073.    not   dx
  1074.    add   ax, 1
  1075.    adc   dx, 0
  1076.    
  1077. CalcExp:
  1078.    div   Ln2Fg16
  1079.    mov   a, ax
  1080.    mov   Num, dx                    ; 0 <= Num < Ln(2)
  1081.       
  1082.    xor   Factorial, Factorial
  1083.    stc
  1084.    rcr   Factorial, 1
  1085.    mov   one, Factorial
  1086.    mov   e, Num
  1087.    mov   Term, Num
  1088.    shr   Num, 1
  1089.       
  1090. Loop_e2x:
  1091.    TaylorTerm
  1092.    add   e, Term                 ; e = x + x*x/2 + (x**3)/3! + . . .
  1093.    cmp   Term, WORD PTR TrigLimit
  1094.    jnbe  SHORT Loop_e2x
  1095.       
  1096. ExitIntSinhCosh:
  1097.    stc
  1098.    rcr   e, 1                    ; e = e + 1, fg 15
  1099.    ret                           ; return e**x * (2**15), 1 < e**x < 2
  1100. _e2x   ENDP
  1101.       
  1102.       
  1103.       
  1104. Exp086    PROC     uses si di, LoNum:WORD, HiNum:WORD
  1105.    mov   ax, LoNum 
  1106.    mov   dx, HiNum
  1107.    
  1108.    call  _e2x
  1109.       
  1110.    cmp   a, 16
  1111.    jae   Overflow
  1112.       
  1113.    cmp   expSign, 0
  1114.    jnz   NegNumber
  1115.       
  1116.    mov   ax, e
  1117.    mov   dx, ax
  1118.    inc   a
  1119.    mov   cx, 16
  1120.    sub   cx, a
  1121.    shr   dx, cl
  1122.    mov   cx, a
  1123.    shl   ax, cl
  1124.    jmp   ExitExp086
  1125.       
  1126. Overflow:
  1127.    xor   ax, ax
  1128.    xor   dx, dx
  1129.    mov   TrigOverflow, 1
  1130.    jmp   ExitExp086
  1131.       
  1132. NegNumber:
  1133.    cmp   e, 8000h
  1134.    jne   DivideE
  1135.       
  1136.    mov   ax, e
  1137.    dec   a
  1138.    jmp   ShiftE
  1139.       
  1140. DivideE:
  1141.    xor   ax, ax
  1142.    mov   dx, ax
  1143.    stc
  1144.    rcr   dx, 1
  1145.    div   e
  1146.       
  1147. ShiftE:
  1148.    xor   dx, dx
  1149.    mov   cx, a
  1150.    shr   ax, cl
  1151.       
  1152. ExitExp086:
  1153.    ret
  1154. Exp086    ENDP
  1155.  
  1156.  
  1157.  
  1158. SinhCosh086    PROC     uses si di, LoNum:WORD, HiNum:WORD, SinhAddr:WORD, \
  1159.                                    CoshAddr:WORD
  1160.    mov   ax, LoNum
  1161.    mov   dx, HiNum
  1162.  
  1163.    call  _e2x
  1164.  
  1165.    cmp   e, 8000h
  1166.    jne   InvertE              ; e > 1
  1167.  
  1168.    mov   dx, 1
  1169.    xor   ax, ax
  1170.    cmp   a, 0
  1171.    jne   Shiftone
  1172.  
  1173.    mov   e, ax
  1174.    mov   cx, ax
  1175.    jmp   ChkSinhSign
  1176.  
  1177. Shiftone:
  1178.    mov   cx, a
  1179.    shl   dx, cl
  1180.    dec   cx
  1181.    shr   e, cl
  1182.    shr   dx, 1
  1183.    shr   e, 1
  1184.    mov   cx, dx
  1185.    sub   ax, e
  1186.    sbb   dx, 0
  1187.    xchg  ax, e
  1188.    xchg  dx, cx
  1189.    jmp   ChkSinhSign
  1190.  
  1191. InvertE:
  1192.    xor   ax, ax               ; calc 1/e
  1193.    mov   dx, 8000h
  1194.    div   e
  1195.  
  1196.    mov   Inve, ax
  1197.  
  1198. ShiftE:
  1199.    mov   cx, a
  1200.    shr   Inve, cl
  1201.    inc   cl
  1202.    mov   dx, e
  1203.    shl   e, cl
  1204.    neg   cl
  1205.    add   cl, 16
  1206.    shr   dx, cl
  1207.    mov   cx, dx               ; cx:e == e**Exp
  1208.  
  1209.    mov   ax, e                ; dx:e == e**Exp
  1210.    add   ax, Inve
  1211.    adc   dx, 0
  1212.    shr   dx, 1
  1213.    rcr   ax, 1                ; cosh(Num) = (e**Exp + 1/e**Exp) / 2
  1214.  
  1215.    sub   e, Inve
  1216.    sbb   cx, 0
  1217.    sar   cx, 1
  1218.    rcr   e, 1
  1219.  
  1220. ChkSinhSign:
  1221.    or    HiNum, 0
  1222.    jns   StoreHyperbolics
  1223.  
  1224.    not   e
  1225.    not   cx
  1226.    add   e, 1
  1227.    adc   cx, 0
  1228.  
  1229. StoreHyperbolics:
  1230.    mov   bx, CoshAddr
  1231.    mov   WORD PTR [bx], ax
  1232.    mov   WORD PTR [bx+2], dx
  1233.  
  1234.    mov   bx, SinhAddr
  1235.    mov   WORD PTR [bx], e
  1236.    mov   WORD PTR [bx+2], cx
  1237.  
  1238.    ret
  1239. SinhCosh086    ENDP
  1240.  
  1241.  
  1242.  
  1243. Numerator   equ      <bx>
  1244. Denominator equ      <di>
  1245. Counter     equ      <si>
  1246.  
  1247. Log086   PROC     uses si di, LoNum:WORD, HiNum:WORD, Fudge:WORD
  1248. LOCAL Exp:WORD, Accum:WORD, LoAns:WORD, HiAns:WORD
  1249.       xor   bx, bx
  1250.       mov   Accum, bx
  1251.       mov   LoAns, bx
  1252.       mov   HiAns, bx
  1253.       mov   cx, Fudge
  1254.       mov   ax, LoNum
  1255.       mov   dx, HiNum
  1256.  
  1257.       or    dx, dx
  1258.       js    Overflow
  1259.       jnz   ShiftUp
  1260.  
  1261.       or    ax, ax
  1262.       jnz   ShiftUp      
  1263.  
  1264.    Overflow:
  1265.       mov   TrigOverflow, 1
  1266.       jmp   ExitLog086
  1267.  
  1268.    ShiftUp:
  1269.       inc   cx                      ; cx = Exp
  1270.       shl   ax, 1
  1271.       rcl   dx, 1
  1272.       or    dx, dx
  1273.       jns   ShiftUp
  1274.  
  1275.       neg   cx
  1276.       add   cx, 31
  1277.       mov   Exp, cx
  1278.       mov   di, dx                  ; dx = x, 1 <= x < 2, fg 15
  1279.       add   di, 8000h
  1280.       rcr   di, 1                   ; di = x + 1, fg 14
  1281.  
  1282.       shl   ax, 1
  1283.       rcl   dx, 1
  1284.       push  ax
  1285.       push  dx
  1286.  
  1287.       mov   ax, dx
  1288.       mul   dx
  1289.       mov   Numerator, dx           ; Numerator = (x - 1)**2, fg 16
  1290.                                     ; 0 <= cx < 1
  1291.  
  1292.       mov   ax, di
  1293.       mul   di
  1294.       mov   si, dx                  ; Save (x + 1)**2, fg 12
  1295.       pop   dx
  1296.       pop   ax
  1297.       shr   dx, 1
  1298.       rcr   ax, 1
  1299.       div   di
  1300.       mov   LoAns, ax
  1301.       push  ax                      ; Rem = (x-1)/(x+1), fg 17
  1302.       mov   Denominator, si         ; Denominator = (x + 1)**2, fg 12
  1303.                                     ; 1 <= Denominator < 8
  1304.       mov   Counter, 3
  1305.       mov   cl, 3
  1306.       mov   ch, 13
  1307.  
  1308.    LogLoop:
  1309.       pop   ax                      ; ax = Rem, fg 13+cl
  1310.       mul   Numerator               ; dx:ax = Rem * (x-1)**2, fg 29+cl
  1311.       shr   dx, 1
  1312.       rcr   ax, 1
  1313.       div   Denominator
  1314.       push  ax                      ; Rem = [Rem*(x-1)**2]/(x+1)**2, fg 17+cl
  1315.       xor   dx, dx
  1316.       div   Counter                 ; ax = Rem / Counter, fg 17+cl
  1317.       mov   dx, ax
  1318.       shr   ax, cl                  ; ax = Rem / Counter, fg 17
  1319.       cmp   ax, WORD PTR TrigLimit
  1320.       jbe   SHORT MultExp
  1321.  
  1322.       xchg  cl, ch
  1323.       shl   dx, cl
  1324.       add   Accum, dx
  1325.       adc   LoAns, ax
  1326.       adc   HiAns, 0
  1327.       inc   Counter
  1328.       inc   Counter
  1329.       xchg  cl, ch
  1330.       add   cl, 3
  1331.       sub   ch, 3
  1332.       jmp   LogLoop
  1333.  
  1334.    MultExp:
  1335.       pop   cx                      ; Discard Rem
  1336.       mov   cx, Exp
  1337.       mov   ax, Ln2Fg16
  1338.       or    cx, cx
  1339.       js    SubFromAns
  1340.  
  1341.       mul   cx                      ; cx = Exp * Lg2Fg16, fg 16
  1342.       add   LoAns, ax
  1343.       adc   HiAns, dx
  1344.       jmp   ExitLog086
  1345.  
  1346.    SubFromAns:
  1347.       neg   cx
  1348.       mul   cx
  1349.       sub   LoAns, ax
  1350.       sbb   HiAns, dx
  1351.  
  1352.    ExitLog086:
  1353.       ; r = (x-1)/(x+1)
  1354.       ; Ans = 2 * (r + r**3/3 + r**5/5 + . . .) + (Exp * Ln2Fg16), fg 16
  1355.       mov   ax, LoAns
  1356.       mov   dx, HiAns
  1357.       ret
  1358. Log086   ENDP
  1359.  
  1360.  
  1361. END
  1362.